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PREDICTION OF AXIAL-FLOW INSTABILITIES IN A TURBOJET ENGINE 


BY USE OF A MULTISTAGE COMPRESSOR SIMULATION 
ON THE DIGITAL COMPUTER 
by Carl J. Daniele, Ronald J. Blaha, and Kurt Seldner 
Lewis Research Center 
SUMMARY 

A method of estimating the undistorted stall line for an axial-flow compressor 
by using the digital computer is presented . The method involves linearization of 
nonlinear dynamic equations (state equations) about an operating point on a speed 
line . The linearization process is accomplished by using finite differences . (Each 
state variable is incremented in turn while all others are held constant . The result- 
ing change in the state equations is calculated .) The changes in the state equations 
are divided by the change in the state variable to give the first partial derivative . 
These first partial derivatives form a linear system matrix . Then the first method 
of Lyapunov is applied to determine the stability of the nonlinear system from the 
stability of the linear system . The method is applied to a simulation of the J 85 com- 
pressor which utilizes stage stacking and lumped volume techniques for the inter- 
stage regions to simulate steady-state and dynamic compressor performance . The 
stability boundary predicted by the digital simulation compares quite well with the 
stall line predicted by a dynamic simulation of the J85 compressor programmed on 
the analog computer . The method of stall prediction on the analog computer is to 
proceed up a compressor speed line until a small amplitude limit cycle is observed. 
Since previous studies have shown that the analog-predicted stall line agrees well 
with the experimental stall line of the compressor , the digital method presented is 
also a good means of estimating the stall line . 



INTRODUCTION 


The steady-state and transient performance of airplane propulsion systems 
using axial-flow compressors is greatly affected by the compressor-stall character- 
istics . Compressor stall adversely affects engine performance and engine life and 
can cause inlet unstarts in supersonic aircraft. Thus, the need for determining the 
stall line for these compressors is evident . 

Until recently stall lines for axial-flow compressors could only be determined 
experimentally . However , engine experiments to define stall lines are costly and 
time consuming . Thus , there is a need to determine these stall lines analytically . 

In reference 1 it is shown that an analog simulation can accurately predict the un- 
distorted stall line for the J85 turbojet engine. However, the ability to predict stall 
depends on the stage characteristics of the compressor that must be determined 
experimentally. The eight- stage axial compressor was modeled on the analog com- 
puter by using stage performance maps and a volume behind each stage to simulate 
gas dynamics (lumped volume) . Also included in the simulation was a burner volume 
followed by a choked nozzle (simulating a choked turbine nozzle) . 

The procedure to define stall on the analog computer is to start at a steady- state 
operating point on a speed line and to slowly decrease nozzle area . This causes an 
increase of back pressure on the compressor so that the operating point moves up 
the speed line . At some point a limit cycle is observed , which increases in amplitude 
as the operating point continues up the speed line . The compressor stall point on 
that speed line is assumed to be the point at which the limit cycle is first observed . 

A digital simulation of the J 85 engine was also formulated . The compressor was 
modeled on the digital computer in the same way it was modeled on the analog . The 
digital simulation transient analysis proved to be too time consumming because of 
the serial and discrete (rather than parallel and continuous) nature of the digital 
computer . To simulate dynamics correctly , the time step used must be smaller than 
the smallest time constant in the system , which necessitated a time step of 50 micro- 
seconds for the J85 simulation. In addition, small amplitude limit cycles are difficult 
to detect and distinguish from numerical instability , which may be induced by the 
numerical integration process ... 

Therefore , a new method for determining the stall line from the digital simula- 
tion was developed and is the subject of this report . The method involves lineariza- 



tion of the nonlinear system of equations (state equations) that defines the compres- 
sor simulation about a steady-state operating point, thereby obtaining a linear sys- 
tem matrix . The elements of the system matrix are the partial derivatives of the 
state equations with respect to the state variables . The linearization of the nonlinear 
equations and the calculation of matrix elements is accomplished by finite differences . 
The eigenvalues of the system matrix are the poles of the system . If the poles all 
have negative real parts (in the left half plane in Laplace space) , the system is 
stable . However , if one or more eigenvalues has a positive real part , the system 
is unstable . The stall point on the speed line is assumed to be the point at which 
positive eigenvalues are encountered . This method for determining the stability of 
the nonlinear system from the stability of the linear system is based on the first 
method of Lyapunov . The procedure is readily mechanized on the digital computer 
and requires very little computer time because there is no integration of dynamic 
variables . 

The report is arranged in the following manner: First a short description of the 
J85 engine is given. Then a complete description of the steady-state and dynamic 
equations used to model the compressor is given . The steady-state equations are 
used to obtain an operating point on a speed line . Next , the methods of implementing 
the compressor model on analog and digital computers and the methods of stall pre- 
diction from each simulation are described . The stall lines predicted by both simu- 
lations are compared . Also compared are the frequency responses of the two meth- 
ods . After this an eigenvector analysis of unstable modes of the linearized system 
is presented . Finally , an approximation is made to simplify the compressor repre- 
sentation and thus reduce the number of state equations . This results in a smaller 
system matrix . The stall line prediction of this simplified model is compared to the 
full variable model . 


ENGINE DESCRIPTION 

The analysis in this report is applied to the J 85-13 afterburning turbojet engine. 
The engine consists of an eight- stage, axial-flow compressor, a two- stage turbine, 
variable inlet and outlet guide vanes, controlled interstage bleeds, and variable 
exhaust nozzle . The inlet guide vanes are scheduled as a function of corrected 
speed; the outlet guide vanes are fixed, but the pressure drop across them varies 
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with corrected weight flow . The interstage bleeds are on the third , fourth , and 
fifth compressor stages and are also scheduled as a function of corrected speed . 

(The bleeds are closed above 95.5 percent corrected speed, are full open at 80 per- 
cent speed , and below , and vary with speed in between . ) A simulation of the whole 
engine is described in reference 2 . 


SIMULATION MODEL 

The multistage compressor portion of the analog and digital simulations uses 
steady-state stage characteristics and an associated stage volume for gas dynamics 
to represent individual stages . Each compressor stage is modeled with pressure 
and temperature maps; each volume, with lumped continuity, energy, and momentum 
equations . The inlet and outlet guide vanes are modeled as pressure drops into and 
out of the first and last stages , respectively - the amount of pressure drop being a 
function of corrected speed on the digital and corrected weight flow on the analog. 
Bleed flows are included for the third, fourth, and fifth stages; the amount of bleed 
flow is also a function of corrected speed . In both analog and digital simulations , 
the bleeds are assumed choked. Between the last compressor stage and turbine, a 
burner volume is modeled . No heat addition is included , however , because only 
the compressor , rather than the whole engine , is of interest for this report , Behind 
the burner a choked nozzle is used to terminate the simulation . The choked nozzle 
simulates a choked turbine in the actual engine . A schematic of an ideal compressor 
stage is shown in figure 1 . All symbols are defined in appendix A , The steady- state 
representation of the compressor stage is presented in appendix B , and the equa- 
tions for the lumped volume gas dynamics are given in appendix C . 


Steady-State Equations 

The steady-state performance of the axial-flow compressor is calculated by 

using the stage stacking technique described in reference 1 . In this technique the 

discharge conditions of one stage are used as the inlet conditions to the next stage . 

The steady- state performance for each stage is given by pressure and temperature 

P T 

rise coefficients v|/ and , These coefficients are defined (in appendix B) as 

C C 
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Representative maps of \|/ and \|/ as a function of flow coefficient <& are shown 
in figures 2 and 3 , respectively . The flow coefficient is defined as 
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The axial-flow velocity v is related to the stage-inlet total pressure and tem- 

ZC 9 XI 

perature (as developed in appendix B) by 

W c> n V e v, n-1 _ v zc, n f / v zc, n 
5 v, n-l^c, n V e v, n-1 \^v, n-1 

where 8 and 0 are corrected pressure and temperature parameters , and p is 
the angle between the air velocity and the axial direction . In steady state the inlet 
and outlet guide vanes are modeled as pressure drops , thus 
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P_ ft » P2 (IGVPR) - 0. 005 P2 (5) 

V, u 

P OGV “ P tv, 8 (OGVPH) (6) 

The steady-state performance of the burner volume is modeled by assuming that 
the pressure , temperature , and weight flow into and out of the volume are the same: 

T tv, 9 = T tc, 9 (7) 

P tv, 9 = p tc, 9 (8) 
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while 


p to,9= P OGV- K '<W c ,9> 2 


(9) 


where K' is determined from experimental data. The nozzle area and bleed flows 
are assumed choked , with the amount of bleed flow area a function of corrected 
speed . The choked flow equation used is 


W = KCHOK s4 


( 10 ) 


The steady- state equations are used to determine an operating point on a speed line. 


Dynamic Equations 

Compressor dynamics are simulated by using a lumped volume for each stage . 
Each lumped volume is located just downstream of the stage’s active element. The 
volume dynamics are modeled through the application of continuity , energy , and 
momentum equations as shown in appendix C . The equations are 

~ (Psv, n) = (We, n “ ^c, n+l) (11) 

at v n 

jj- (Psv, n T tv, n) = [ T tc, n W C) n “ Ttv, n (W c , n+1 + ^b, n)] (12) 

and 

* n) = ( p to, n - p tv, n) < 13 ) 


While the equation of state for the gas in the volume is 

P tv, n = R Psv, n T tv, n (14) 

No dynamic effects are modeled for the inlet and outlet guide vanes , interstage 
bleeds , or choked nozzle . A lumped volume is included in the burner but without 
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an energy equation. Instead, equality of temperature is assumed: 


T tv, 9 T te, 9 T tv, 8 


( 15 ) 


Analog Implementation of Equations 

Equations (1) to (15) form a complete mathematical representation of the dynam- 
ic and steady- state characteristics of the compressor model. These equations are 
programmed on the analog computer. All functional relationships for steady-state 
stage characteristic pressures and temperatures and for bleed flows and guide 
vanes are generated by using diode function generators . Dynamics on the analog 
computer are calculated by integrating the dynamic equations (11) to (13) . The 
temperatures in the stage volumes are calculated by integrating equations (11) 
and (12) to obtain p and pT, respectively. 

The method of stall prediction on the analog simulation is to start with a steady- 
state operating point on a speed line . The nozzle area is then slowly decreased so 
that back pressure on the compressor is slowly increased . This causes the operating 
point to move slowly up the speed line . At some point a small amplitude limit cycle 
is observed . Further movement up the speed line results in an increase in amplitude 
of the limit cycle . The stall point for the speed line is assumed to be the point at 
which the limit cycle is first observed . By repeating this procedure for a number 
of speed lines , the undistorted stall line is defined . 


Digital Implementation of Equations 

The equations could also be programmed into the digital computer . It would 
be possible to determine the stall line for the compressor on the digital computer 
in the same way it was done on the analog simulation . Digital simulation of a dynam- 
ic system with small time constants necessitates a time step smaller than the smallest 
time constant in the system . This results in large amounts of computer time to obtain 
a transient solution . In addition , small amplitude limit cycles generated on a digital 
simulation are difficult to detect and distinguish from numerical instability induced 
by the integration method . The same steady-state and dynamic equations used on 
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the analog computer are also used on the digital computer to describe the system . 
However, on the digital computer, the dynamic equations are not integrated, and 
no transient solution is obtained . 

STABILITY ANALYSIS 

Consider the system of nonlinear differential equations: 


X = f(X, u*) (16) 

where X and f are Z -vectors (Z equals 26 for this problem) and u* is an m- 
vector of inputs . One way to determine the stability of this nonlinear system is the 
method of Lyapunov , as presented in reference 3 . Let the steady-state solution of 
this system of equations be given by X = X Q , so that f(X Q ,uJ) = 0. This system can 
be linearized and the steady- state operating point translated to the origin; the 
resulting equations are given by 

Y = A Y + Bu (17) 


where 


Y = X - X, 


u = u* - u* 
o 


A is the Jacobian matrix , such that 

A t'-hi 

ij “ 8Xj j = 1, l 

and B is an l by m matrix, such that 

i=l, l 

B ij = 3Uj j = 1, m 


(18) 


(19) 


( 20 ) 
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For the free dynamic system Y = AY , Lyapunov showed that , if the eigenvalues 
X of A have nonzero real parts , then the stability of the operating point X for 
the nonlinear system is the same as that of the operating point Y = 0 for the linear- 
izied system . The eigenvalues are the roots of the characteristic equation for the 
linearizied system and are determined by 

| A - Al| «= 0 (21) 

If all the eigenvalues have negative real parts (lie in the left half plane in Laplace 
space) , the linear system is stable , and hence the nonlinear system is stable in a 
region about the operating point . If at least one eigenvalue has a positive real part , 
the nonlinear system is unstable. 

The first method of Lyapunov is applied to the compressor simulation already 
described . An operating point is chosen on a speed line and the steady-state values 
for the state variables are determined . The 26 state equations are linearized about 
the operating point , and the system matrix A of first partial derivatives is gener - 
ated . This is accomplished by using finite differences . Each state variable is varied 
in turn , while the others are held constant . If each state variable is changed by d.. , 
then: 
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W x j 


0 (zero in steady state) 


( 22 ) 


The eigenvalues of A are calculated by use of digital programs presented in ref- 
erence 4 . If all the eigenvalues have negative real parts , the system is stable at that 
operating point . (The compressor is not in stall . ) A point is chosen farther up the 
speed line, and the same procedure is performed. This process is continued until 
at least one eigenvalue has a positive real part (in the right half plane in Laplace 
space) . At this point the compressor is assumed to be in stall, and the stall point 
on the speed line is assumed to be the point just before the positive eigenvalue is 
encountered . 

A comparison of analog and digital stability assumptions is of interest . For 
the analog simulation stall is assumed to occur when a limit cycle is observed . For 
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the digital simulation stall is assumed to occur when positive eigenvalues are en- 
countered. It can be shown rather easily and is implied by the Lyapunov theory 
that, if the linear system is stable, the states converge to their steady-state values 
for the nonlinear system. Thus, if all the eigenvalues of the linear system have neg- 
ative real parts, there can be no limit cycle for the nonlinear system. However, if 
the linear system is unstable , the states of the nonlinear system may diverge or may 
have a limit cycle . 

RESULTS 

Stall Line 

The comparison of stall points generated by the analog and digital J85 compres- 
sor simulations is shown in figure 4. Ten corrected speed lines are shown. The 
speed lines shown were obtained on the analog computer . (Speed lines obtained on 
the digital computer were nearly identical.) The analog generated and digital gen- 
erated stall points shown in figure 4 are quite close . The slight error in the 85 and 
87 percent speed lines can be attributed to the difficulty in defining the exact point 
where a limit cycle started with the analog simulation . The results in figure 4 are 
for sea-level- static inlet conditions. However, since no Reynolds number effects 
are included in the simulation , the results apply only for high Reynolds number 
conditions . The results show that the digital-state variable analysis can duplicate 
analog results. The analog method of determining the undistorted stall line has al- 
ready been verified . Thus , the digital state variable analysis is also a good means 
of determining the undistorted stall line for axial-flow compressors . 

Figure 5 shows a pseudo-root locus for the digital compressor representation. 

It shows the motion of eigenvalues through Laplace space as the operating point ap- 
proached stall on the 100 -percent speed line. This is not a conventional root locus 
since the whole system is essentially redefined at each point as the operating point 
moves up the speed line. The figure shows only the upper half plane since the 
lower half is the mirror image . The log scale is used since large decade changes 
occurred in the eigenvalues that did move significantly. However, since there is 
no zero on a log scale , the eigenvalues must be assumed to go through zero some- 
where between the values -0.001 and +0.001 on the abscissa. The circle data 
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symbols indicate that very little change in the eigenvalue position occurred as the 
operating point changed . The five operating points shown started just before stall 
(points 1 to 3) where point 3 is the stall point given in figure 4 . Operating point 4 is 
just into stall and operating point 5 is farther into stall . The decrease in weight flow 
from one operating point to the next was 0.045 kilogram per second. Note that for 
the speed line shown four eigenvalues (two shown plus their complex pairs) crossed 
the imaginary axis almost simultaneously . To determine which pair crossed first , a 
smaller weight flow step would have to be taken; however , this was not done because 
a tolerance of 0 . 045 kilogram per second is within experimental accuracy . The fifth 
operating point is shown because it indicates that further movement up the speed 
line caused those eigenvalues already in the right half plane to move farther to the 
right and another eigenvalue to move closer to becoming unstable . 


Frequency Response 


The frequency responses generated by the analog and digital simulations were 
also compared. The two inputs available for the frequency response were P2, the 
input pressure and , the nozzle area . The analog response was measured by 
using a discrete frequency analyzer with the simulation . The digital response was 
calculated by taking the Laplace transform of equations (17) . Matrix B was calcu- 
lated by keeping all the state variables constant and perturbing each input separately 
while keeping the other constant. The resulting change in Y was calculated and 
the B matrix determined as 


AY 

AY 

1 

1 

Au^ 

Au 2 


AY 

2 

Au-j^ 

Au 2 


\Y 

26 

A *26 

Au.^ 

< 


(23) 
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This procedure was accomplished on the digital computer by using finite differences . 
Then , using matrices A and B the frequency response was calculated using the digi- 
tal programs found in reference 5 . A comparison of frequency responses of the analog 
and digital compressor simulations is shown in figures 6. The operating point chosen 
for the frequency responses was on the 100-percent speed line at a weight flow of 
19.739 kilograms per second and a pressure ratio of 6.79. Three frequency re- 
sponses are shown - P v ^/P2 , W c ^/P2 , and p qGV^ 10 ' In 311 cases there was good 
agreement up to about 20 hertz, at which point the two responses started to diverge. 
The reason for the discrepancy, however, is not clear since both simulations con- 
tained the same volume dynamics . 

Eigenvector Analysis 

Digital programs from reference 4 also have the capability of determining the 
eigenvectors associated with the eigenvalues . The eigenvectors are defined by 

Ae i = Ae i (24) 

where eigenvector e. is a nonzero vector and is associated with eigenvalue X . . 

1 X.t 1 

From reference 3 , if X. is real, the motion e 1 is along the eigenvector e. . If 

X.t 1 

is greater than zero , then e 1 along e. approaches infinity; if less than zero, 

it approaches the origin of the state . If X. is complex , the combined motion of it 

and its complex conjugate is oscillatory and takes place in the plane determined by 

their associated eigenvectors . If the real part of the complex eigenvalue is negative, 

the oscillation is damped; if positive, the oscillation goes unstable. 

The eigenvectors were calculated for the operating point just beyond stall, 
where at least one eigenvalue had a positive real part. For each eigenvalue there was 
a corresponding eigenvector (26 for this system) . Each eigenvector had 26 com- 
ponents - the first component being related to the first state variable and so on. 

Since the eigenvectors represent the fundamental modes of the system , the state 
variables p and pT were changed to pressure and temperature, using equations 
(11) , (12) , and (14) , to relate the results to more physically understandable quan- 
tities . For all 10 speed lines , the eigenvalues that were unstable were complex; thus, 
there were two eigenvectors associated with the instability . The magnitudes of the 
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components of the unstable eigenvectors were calculated and divided by the magni- 
tude of the largest component of each quantity, respectively. (All pressures changes 
were divided by the largest pressure change; all weight flow changes by the largest 
weight flow change; and all temperature changes by the largest temperature change . ) 
This nondimensionalized the numbers so that they could be compared . The resulting 
modal maps for the 100-percent speed line are shown in figure 7. For this speed 
line the values given are for the unstable eigenvalue at an imaginary value of 
11 000 radians per second in figure 5. Note from figure 7 that the largest relative 
motion in pressure occurred in the eighth stage , weight flow in the third stage , 
and temperature in the second stage . Although the phase is not shown , pressure and 
weight flow motions are 90° out of phase with weight flow leading pressure . Tem- 
perature motion is in phase with pressure. Also, for any two successive stages, 
vectors of like quantities could be 180° out of phase. It is impossible to state which 
stage caused stall since the eigenvectors are really a composite of the total motion 
in the system . 

No conclusions could be drawn from these results , as the significance of this 
method has not been fully investigated . Different state variables or a closer look 
at actual stage pressure and temperature ratios may give more definitive results . 

Modified Simulation 


The digital simulation consists of three state variables per stage . However , 
if it is desired to predict the loss in performance and surge margin due to distorted 
inlet flow , circumferential or radial effects would have to be included in addition 
to axial flow , and the number of state variables could increase drastically . Alterna- 
tively , it may be desired to determine the undistorted stall line for a more complex 
engine (such as a two spool turbofan) , which would also cause a large increase in 
the number of state equations. In anticipation of such problems, a simplification 
of the volume dynamics was made by arbitrarily eliminating the temperature dy- 
namics in the compressor stages: 

^tv, n = T tc , n (25) 

Using equation (25) and combining equations (12) and (14) result in 
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dt <I?tT - n) ■ 


( 26 ) 


r BT tv. n 


(W c> n - W 0> n+1 - W b „) 


n 


I- (We, n) - (Ptc, n - Ptv, n) (13) 

From equations (13) , (25) , and (26) , the modified digital simulation consisted of 
18 state variables - two for each compressor stage and two for the burner volume . 
Here also a choked nozzle was used to terminate the burner volume . 

Figure 8 shows the comparison of stall points generated by the unmodified 
digital and modified digital simulations . The same speed lines are shown as in fig- 
ure 4 . From figure 8 the results compare quite well for the lower speed lines (80 per- 
cent and below) , but the results diverge at the higher speeds. This would indicate 
that the temperature dynamics tend to stabilize the system , the effect being most 
pronounced at the higher speeds . For comparison the analog simulation was also 
modified to eliminate the temperature dynamics . Figure 9 represents a comparison 
of the stall points for the same speed lines using both modified simulations . The 
results show good comparison between the two simulations except at 70 percent 
speed . 

One final method of looking at the effect of the temperature dynamics is shown 
in figure 10 . In this figure is a comparison of the eigenvalue position for the modified 
and unmodified digital compressor simulations just after the stall limit for the modi- 
fied simulation . At the same operating point on the 100-percent speed line , all the 
eigenvalues for the unmodified simulation are to the left of those for the modified 
simulation . Because this figure is in Laplace space , a more stable system for the 
unmodified simulation is indicated . From figure 10 each modified eigenvalue corre- 
sponds to an eigenvalue from the unmodified system . The modified eigenvalue , which 
is unstable, has an imaginary value of 36 000 radians per second. However, from 
figure 5 this eigenvalue remains in about the same relative position as the operating 
point moves up the speed line and does not become unstable . Instead , the two eigen- 
values with imaginary parts of about 11 000 and 15 000 radians per second are the 
ones that become unstable for the unmodified simulation, while remaining stable 

for the modified simulation . Finally , the modified state variable method does give 
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good stall-line prediction for lower speeds , but eliminating the temperature dynam- 
ics does degrade stall prediction at the higher speeds. Therefore, the use of this 
method in solving future problems should be considered carefully , 


CONCLUSIONS 

A new method of determining the undistorted stall line for an axial-flow com- 
pressor simulation using the digital computer has been presented . The simulation 
uses the stage stacking and lumped volume techniques first used to define the stall 
line by the analog computer . The method employs small changes about an operating 
point to form linearized state equations for the system . The stability of the nonlin- 
ear system is obtained from the stability of the linear system by the first method of 
Lyapunov . 

This procedure was applied to an eight-stage simulation of the J85 compressor. 
Three state variables were assumed for each stage . The analysis was done for 10 
compressor speed lines , and the stall point was defined on each . There was good 
agreement between the analog and digital stall-point predictions . Since the analog 
previously matched experimental results , the digital method presented in this report 
is also a good means of determining the stall lines for axial-flow compressors . Also, 
frequency responses from the two simulations were compared . Good agreement be- 
tween the two existed up to about 20 hertz, at which point they started to diverge. 
Using the eigenvalues at a point just into stall, the corresponding eigenvectors 
were calculated . From the eigenvectors , modal maps of pressure , weight flow , and 
temperature were drawn. The largest relative motion in pressure occurred in the 
eighth stage; weight flow, in the third stage; and temperature, in the second stage. 
However, no conclusions could be drawn from this analysis, and further investiga- 
tion is warranted . 

Next , the number of state variables for the system was reduced by arbitrarily 
eliminating the temperature dynamics . This modified method showed good corre- 
lation with the analog for stall prediction at low speeds (80 percent and below) , 
However , the stall point predictions diverged at the high speeds , Evidently , the 
temperature dynamics tended to stabilize the system , and this stabilizing effect was 


15 


4 



more pronounced at the higher speeds . Thus , reducing the number of state variables 
in future work should be considered carefully . 

Lewis Research Center, 

National Aeronautics and Space Administration , 

Cleveland, Ohio, July 25, 1974, 

501-24. 


16 



APPENDIX A 


SYMBOLS 

A l by 7 Jacobian matrix 

2 2 

si area, m ; ft 

a speed of sound, m/sec; ft/ sec 

B l by m input coefficient matrix 

C coefficient, rpm 2 /K 2 ; rpm 2 /°R 2 

C specific heat at constant pressure , J/ (kg-K); Btu / (lbm-°R) 

r 

C v specific heat at constant volume, J/ (kg-K); Btu/(lbm-°R) 

d perturbation of state variable 

e eigenvector 

2 2 
g gravitational constant , 1.0 (kg-m)/ (N-sec ); 32.17 (lbm-ft)/ (lbf-sec ) 

h enthalpy, J/kg; Btu/lbm 

I identity matrix 

IGVPR pressure drop due to inlet guide vanes 

J mechanical equivalent of heat , 1.0 (N-m/J); 778.3 (ft-lbf/Btu) 

K coefficient, rpm-sec/m; rpm-sec/ft 

K' constant for burner volume , (N/m 2 )/ (kg/ sec) 2 ; Obf/ft 2 )/ (lbm/sec) 2 

KCHOK choked flow coefficient , (kg-^Kj/ (N-sec); (lbm-^Itj/ (lbf-sec) 
k bleed flow coefficient , |kg-^C^/ (N-sec); ^lbm-V°R^/ (lbf-sec) 

L length, m; ft 

l total number of states 

M Mach number 

N rotational speed, rpm 

OGVPR pressure drop due to outlet guide vanes 
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pressure, N/m 2 ; lbf/ft 2 

. 2 9 

input pressure , N/m ; lbf/ft 

internal energy, J/kg; Btu/lbm 

universal gas constant , 287 (N-m)/ (kg-K); 53.3 (lbf-ft)/ (lbm-°R) 
radius, m; ft 
temperature , K; °R 
time, sec 

rotor speed, m/sec; ft/sec 

m by 1 vector of system inputs translated to the operating point 

m by 1 vector of system inputs 

volume, m 3 ; ft 3 

velocity, m/sec; ft/sec 

weight flow , kg/sec; lbm/sec 

l by 1 vector of state variables 

l by 1 vector of state variables translated to operating point 
coefficient, 30 (m-sec)/ (m-min); 360 (in.-sec)/(ft-min) 
rotor air inlet angle , deg 
ratio of specific heats 
change in variable 

ratio of total pressure to standard atmospheric pressure 

ratio of total temperature to standard atmospheric temperature 

eigenvalue, rad /sec 

weight density, kg/m 3 ; lbm/ft 3 

flow coefficient 

pressure coefficient 

temperature coefficient 



Subscripts: 

b variable associated with stage bleed 

c variable associated with stage characteristic 

i matrix component position 

j matrix component position 

n stage number designation 

OGV outlet guide vanes 

o steady- state operating point 

r reference state 

s static condition 

t total condition 

v variable associated with stage volume 

z axial direction 

0 input to stage one 

1 , 2 , . . . , 8 stage numbers 

9 burner volume 

10 nozzle exit 
Superscripts: 

derivative with respect to time 
— maximum value of quantity used for normalizing 



APPENDIX B 


STEADY-STATE COMPRESSOR REPRESENTATION 


The analytical steady-state model of the compressor is based on the stage 
stacking technique , The derivation of the equations is taken from appendix B of 
reference 1 and is presented for completeness . The stage characteristics are rep- 
resented by a pair of two-dimensional curves for each stage. The individual stage 

P 

performance, described by a pressure coefficient \|r and a temperature coefficient 
T 

v|/ , is assumed to be a function of the flow coefficient O . The coefficients are math- 
ematically described by 

p 2gJ AhQ 


2gJ Ah t 

r * — — 

IT 
v 


$ 


_z 

U 


J 


(Bl) 


xi_ 

The ideal total enthalpy rise for the n n stage can be related to the stage total- 
pressure ratio from 


Ah tc, n “ C p T tv, n-1 


( 7 - l)./y 


tc, n 
tv, n-1/ 


- 1 


The adiabatic enthalpy rise can be determined from 





(B2) 


(B3) 
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and the stage wheel speed U from 


U = 


7rNr 

oi 


(B4) 


The stage inlet conditions are associated with the conditions in the preceding vol- 
ume, and the variables are subscripted accordingly . The stage exit conditions are, 
however, associated with the stage itself. Assuming that the n n compressor stage 
is preceded by the n - 1 volume , we can combine equations (Bl) with equations 

xu 

(B2) to (B4) to form the following set of equations for the n n stage (assuming 
y = 1.4): 


P s £n_ T tv, n - l 

N 2 


2/7 


'\ 


tc, n 
^ p tv, n-l) 


- 1 


T _ C n 

“ n 2 < T tc, n - T tv, n) 




(B5) 


K n 

#n ~ Y V ze, n 


where 




2gJC p (q) : 

n 2 r 2 

n 


K «= Q? 
n-7rr n 


The flow coefficient is defined as the ratio of the axial-flow velocity to a mean 
rotor velocity . Introducing the square root of the stage corrected temperature into 
the relation for <t> given in equations (Bl) , the flow coefficient becomes 
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vV«~ 

$ = _ 

u/vr 

The axial velocity v can be computed from the flow , total pressure , and total 

z 

temperature at the stage inlet . The flow and the axial velocity are related by 


z 


Eliminating the density from equation (B7) with a state equation of the form 

P s = Ps R T s 


results in 


_ WR T s 
T *° * p s 




1 /(^ 1 ) 


Equation (B9) can be rewritten in terms of corrected parameters as 




s 1/(7- 1) 




•1/(7- 1) 
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It can be shown from the energy equation that 


or 



££p 

yR 





Thus, since 



combining equations (BIO) and (B12) results in 


wy£ v z 

s/8 “ye 



1 /( 7 - 1 ) 

Ptr 


It should be noted that the nonsub scripted velocity v introduced into equation 
is the magnitude of the velocity vector and not the axial velocity component v 

z 

If the airflow approaches the rotor inlet angle |3 , 

v z = v cos /3 


and equation (B13) becomes 


wVe = _^z 

s/8 yj 


i 



1/(7- 1) 


1 

2gJC p T tr cos 2 /3 


Ptr 


(BH) 


(B12) 

(B13) 

(Bll) 

(B14) 

(B15) 
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and the corrected axial velocity 
of equations (B15) rewritten for 



where 


24 



APPENDIX C 


STAGE GAS DYNAMICS 

The gas dynamics of each compressor stage is assumed to occur in an equivalent 
stage volume , which is conceptually located just downstream of the stage pressure 
and temperature -rise characteristics. The stage gas dynamics are represented by 
continuity , energy , and momentum balances across the volume . The equations for 
the gas dynamics are developed in appendix C of reference 1 and are presented here 
for completeness . Assuming no work is performed in the volume and that the net 

+Vi 

heat input is negligible , an energy balance can be written for the n volume as 

Yn ~ (Psv, n^tv, n) = ^tc, n Wc, n - btv, n ^c, n+1 (Cl) 

Noting that q t = C v T^ and that h t = C p T t > equation (Cl) becomes 

~[j~ (Psv, n T tv, n) = ( T tc, n ^c, n ” p tv, n w c, n+l) ( C2 ) 

n 

To determine the density of the gas in the stage volume , it is assumed that the rate 
of change of the static density is proportional to the net weight flow: 

£ (Psv, n> = (Wc, n - < n+l> <C3> 

n 

The total pressure in the stage volume is calculated from the equation of state of the 
form 

5/2 

p tv,n* R ( 1+0 - 2M ^n) (Psv,n T tv,n> (C4) 

A momentum balance written across the stage and its associated volume is used 

to provide an expression for the determination of weight flow . It is assumed that 

flow acceleration is caused by the difference between the static compressor head 

rise AP and the head rise (P -P ,) imposed by the adjacent volumes . 

sc,n sv,n sv,n-l J 
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Neglecting convective momentum terms , this assumption yields 




dt 


(Wo, n) ~ [^ P sc, n “ ( p sv, n “ p sv, n-l)j 


(C5) 


The stage inlet conditions , however , are conditions in the inlet volume and , thus , 
AP se,n = P sc,n ‘ P sv,n-X and equation < C5 > beeomes 


dt 


(W c , n) 


Ln 


( p sc, n " p sv, n) 


(C6) 


Assuming the Mach number at the stage exit and the volume are the same , equation 
(C6) can be written as 


d • ^ n g 

— (W c rj) * — 
dt v c>iy L n 


1 + 0. 2 M, 


2 \' 7/2 

[ v, nj ( p tc 5 n “ p tv, 


(C7) 


Equations (C2) to (C4) and (C7) provide a representation for the gas dynamics 
of the n in stage volume . In the computer simulation of these equations (both analog 
and digital simulations) the distinction between static and total pressures was 
dropped . This simplification resulted in neglecting the Mach number corrections 
in equations (C4) and (C7) , The interstage Mach numbers encountered in the com- 
pressor are less than 0.5. The range of corrections in equations (C4) and (C7) are 
thus 


1.0 < (l+ 0.2 M 2 ) < 1. 13 


5/2 


0. 84 < (l + 0. 2 M 2 ) 7/ " 2 < 1. 0 


The worst-case error introduced by this simplification is , thus , approximately 
equivalent to a 15-percent change in the coefficients . With the computer formulation 
used , the approximation is reflected in dynamic coefficients and does not alter 
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steady-state results . Neglecting the Mach number corrections results in 


V~ < W e, n - W c> n+1 ) * 0 
n 


r 


v. 


n 


(T tc> n W c , n - Ttv, n W c , n+l> ' ® 


and 




L (Ptc, n “ Ptv, n) ~ ® 


in the steady state. Equations (C8) thus yield 


w e, n = W c , n+l 


Ttc, n ~ ^tv, n f 
Ptc, n = Ptv, n 


(C9) 


in the steady state . Thus , the steady-state results are unaffected by approximations 
in the volume dynamics . 

Interstage bleed can be included in the compressor simulation by modifying 
equations (C2) and (C3) . The bleed flow can be computed from the relation: 


Wh n-lctX n - V' p 
V T tv, n 


(CIO) 


which assumes choked flow . Including the effects of interstage bleed changes equa- 
tions (C2) and (C3) to 


d_ 

dt 


(Psv, n T tv, n) = [Ttc, n W c , n “ T tv,n (W c , n+l + Wb, n)J (Cll) 


(Psv, n) “ ■“ (W c , n ~ W c , n+l ~ w b, n) 


dt 


(Cl 2) 


n 
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Figure 6. - Comparison of frequency responses generated by analog and digital J85 
simulations. 
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Figure 9. - Comparison of stall points generated ana- 
lytically by modified analog and digital compressor 
simulations. 1 Sea level static conditions. 
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Figure 10. - Comparison of eigenvalues for modified and unmodified compressor digital simulations 
at point just after stall limit for modified compressor simulation. 
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